This writing is a complimentary document to a power session workshop that was taught by Richard Johansen and Mark Chalmers during Data Day 2019 at the University of Cincinnati. It is written with the expectation that the user is familiar with using Rstudio. If you follow along with these instructions, you will be able to create better, interactive maps of social vulnerability at the county level. For a more in depth explanation as to how the data was retrieved, cleaned, and manipulated, please refer to the full R script called Mapping_Social_Vulnerability.R located in the Scripts folder of the GitHub repository. At this point, you might be asking yourself “What GitHub repository?”.
The original script, all of the data, the powerpoint presentation, and supplemental resources used in the workshop can be downloaded from the GitHub repository located at: https://github.com/RAJohansen/DataDay2019
The first step towards creating these maps yourself is to clone that repository, extract the files, and open the R project file titled DataDay2019. This should open with Rstudio automatically. If it does not know to open with Rstudio, you will have to tell it to use that program.
Using the project file allows you to use the exact same code that is in this script to read in the raw or clean data files that were used in the workshop and the following examples.
The next step is to make sure the following packages are downloaded and installed in your R environment. To install these packages, run the following line of code. You might already have some (such as dplyr and tidyverse) and you should modify this to only download the packages you do not already have. We will not necessarily be using all of these packages in the code included here, but if you want to follow each step in the original script you will need these packages. If you only want to produce maps with the data that is provided here, you will need the following packages: sf, tmap, and leaflet.
install.packages(c("tigris","tmap","tidyverse","tabulizer","dplyr","sf","leaflet"))
After the packages are installed, the next step is to use the library function to make sure we can use the functions these packages provide us. Run the following lines of code and we are ready to get started!
library(tigris)
library(tmap)
library(tidyverse)
library(tabulizer)
library(dplyr)
library(sf)
library(leaflet)
We extracted the data directly from the PDF on the website, located at “http://artsandsciences.sc.edu/geog/hvri/sites/sc.edu.geog.hvri/files/attachments/SoVI_10_14_Website.pdf”.
The following code comes from the Mapping_Social_Vulnerability.R script and shows how to use the extract_tables function to scrape the data from the pdf at that website. This is included as an example, but moving forward we will just be reading in the data without going into the retrieval and cleaning process. If you are curious about that process, refer to the full Mapping_Social_Vulnerability.R script. You do not need to run the code in the chunk below.
Sovi_table <- extract_tables("http://artsandsciences.sc.edu/geog/hvri/sites/sc.edu.geog.hvri/files/attachments/SoVI_10_14_Website.pdf")
# Lets use two more functions to convert the extracted table
# into a more usable and analysis friendly format
final <- do.call(rbind, Sovi_table[-length(Sovi_table)])
# Reformate table headers by dropping the first row
final <- as.data.frame(final[2:nrow(final), ])
# Lets lable the column names so they can merged with Census data
headers <- c('GEOID', 'State_FIP', 'County_FIP', 'County_Name', 'CNTY_SoVI',
'Percentile')
# Apply our names to the data frame
names(final) <- headers
# **NOTE** GEOID is the ID code for CENSUS data
# This is mandatory for the next section
### Step 4: Save the table as a csv
# This is helpful for eliminating redundancy and reproducibility
write.csv(final, file='Data/SoVI.csv', row.names=FALSE)
This data is already saved in the repository you downloaded. Instead, just run this line of code and you will read in the clean data file into your environment.
df <- read.csv('Data/SoVI.csv')
The first five rows of the data frame should look like this:
## GEOID State_FIP County_FIP County_Name CNTY_SoVI Percentile
## 1 1005 1 5 Barbour -0.93 34.8
## 2 1027 1 27 Clay 0.15 52.8
## 3 1091 1 91 Marengo 1.40 73.2
## 4 1069 1 69 Houston 0.19 53.4
## 5 1019 1 19 Cherokee -1.46 26.9
## 6 1003 1 3 Baldwin -1.03 33.5
Now that we have the Social Vulnerability data, lets move forward and make some maps. Instead of looking at County level information for the entire country, we are going to narrow the scope down to one state. Our example state is Florida but the user could look at any state they please by making very slight modifications to the code.
The first line of code reads in the census data from Florida. The second line takes only the social vulnerability data from our entire data set with State_FIP = 12, which translates to Florida. The next line merges these two data sets together.
Counties_FL <- st_read('Data/Counties_FL.gpkg')
df_FL <- subset(df,State_FIP == "12")
FL_SoVI <- merge(Counties_FL,df_FL, by = "GEOID", all = FALSE)
This line of code uses base R’s most basic plotting function to make our first map.
plot(FL_SoVI[1])
Now we are going to use the tmap package, which uses the same grammar of graphics as ggplot. The tmap package is specially designed for visualizing geographic and spatial data.
The first map will be basic and then we will add layers.
tm_shape(FL_SoVI) +
tm_fill()
Let’s add county borders.
tm_shape(FL_SoVI) +
tm_borders() +
tm_fill()
Lets add our Social Vulnerability data to explore trends.
tm_shape(FL_SoVI) +
tm_borders() +
tm_fill(col = "CNTY_SoVI")
Manually define lable breaks.
breaks = c(-6,-3,0,3,6)
tm_shape(FL_SoVI) +
tm_borders() +
tm_fill(col = "CNTY_SoVI",breaks = breaks)
In this case, negative values correspond to higher resilience, or less vulnerability. Lets flip the color scheme to better match our intuition.
tm_shape(FL_SoVI) +
tm_borders() +
tm_fill(col = "CNTY_SoVI",breaks = breaks, palette = "-RdYlGn")
Lets choose our own color palette and add a continuous scale bar
tm_shape(FL_SoVI) +
tm_borders() +
tm_fill(col = "CNTY_SoVI", style = "cont", palette = "viridis")
Add some cartographic elements
tm_shape(FL_SoVI) +
tm_borders() +
tm_fill(col = "CNTY_SoVI", style = "cont", palette = "viridis") +
tm_layout(title = "Florida SoVI Vulnerability Index by County",
legend.outside = FALSE,
frame = TRUE,
inner.margins = 0.1,
legend.title.size = 1.5,
legend.text.size = 1.1) +
tm_compass(type = "arrow", position = c("right", "top"), size = 2) +
tm_scale_bar(breaks = c(0, 100, 200),size = 0.8)
The tmap package allows us to make some much more sophisticated visualizations than base plot. Now that we now how to create these beautiful maps, lets add the extra layer of interactivity to them so that users can better engage with the visualization.
First, convert the static tmap to an interactive map using the leaflet package. The popup.vars argument in the tm_fill function allows us to define what information pops up when the user clicks on a given county. Go ahead and click on a county in the map below, it should give you the specific Social Vulnerability value for that county.
map <- tm_shape(FL_SoVI) +
tm_borders() +
tm_fill(col = "CNTY_SoVI",
palette = "-RdYlGn",
id = "NAME",
popup.vars = c("NAME","CNTY_SoVI"))
tmap_leaflet(map)
Next we are going to create a more complex map with leaflet. First, lets read in the Florida population data that is included in the repository and merge it with the spatial data.
FL_pop <- read.csv("Data/FL_Population.csv")
FL_pop$NAMELSAD <- FL_pop$County
FL_SoVI <- merge(FL_SoVI,FL_pop, by = "NAMELSAD", all = FALSE)
We also included some very coarse sea level predictions in the repository. Let’s read these into the working environment so that they can be added to the visualizations.
FL_slr_10ft <- st_read("Data/FL_slr_10ft.gpkg")
The following code creates duplicate maps so we can do a side by side comparison of Social Vulnerability and population with the projected sea level rise overlayed.
facets <- c("CNTY_SoVI","Population")
map_facets <- tm_shape(FL_SoVI) +
tm_borders() +
tm_fill(col = facets,
palette = "-RdYlGn",
id = "NAME",
popup.vars = c("NAME","CNTY_SoVI", "Population")) +
tm_shape(FL_slr_10ft) +
tm_polygons(col = "blue", alpha = 0.75) +
tm_facets(nrow = 1, sync = TRUE, free.scales.fill =TRUE)
tmap_leaflet(map_facets)
These code shows you how to change the basemap that is underneath the interactive map.
map_facets_base <- tm_basemap(leaflet::providers$Esri.WorldImagery) +
tm_shape(FL_SoVI) +
tm_polygons(facets) +
tm_borders() +
tm_fill(col = facets,
id = "NAME",
palette = "-RdYlGn",
popup.vars = c("NAME","CNTY_SoVI", "Population")) +
tm_shape(FL_slr_10ft) +
tm_polygons(col = "blue", alpha = 0.75) +
tm_facets(nrow = 1, sync = TRUE)
tmap_leaflet(map_facets_base)